//Industry
args ind

/*****************************************************************
Formdata
*****************************************************************/
/*
This program simply applies the code given to me by Salvador Navarro
to form the Chilean data.

The only argument is the sector we want to keep
*/


cap program drop formdata
program formdata
		args ind
		
		
		keep if ciiu_3d==`ind'

		gen adv=.
		replace adv=0 if adverts==0 & adverts~=.
		replace adv=1 if adverts>0 & adverts~=.

		gen skwages=wageswc+ bonuswc+ prtaxwc+ fataxwc
		replace skwages=. if skwages<=0
		gen unskwages=wagesbc+ bonusbc+ prtaxbc+ fataxbc
		replace unskwages=. if unskwages<=0
		gen totwages=skwages+unskwages
		replace totwages=. if totwages<=0
		gen ratiowcbc = whitelab / bluelab
		replace ratiowcbc = 10000 if ratiowcbc == . & bluelab == 0 & whitelab ~= 0

		gen output= salegds + shipgds + resales - goodspch + contract + conrepai + coninsta + prodbld + prodmach + prodveh + elecsval + finvufgd + finvfgds - iinvufgd - iinvfgds
		replace output=. if output<=0
		gen matls = totrawma - finvrm + iinvrm + watival + watdval + greasval + lubeval
		replace matls=. if matls<=0
		replace realmats=. if realmats<=0
		gen energ = coalval + carbval +cokeval + petrval + diesval + benzval + parafval + lgasval + pgasval + fwoodval + ofuelval + elecbval
		replace energ=. if energ<=0 
		replace renerg=. if renerg<=0 
		gen serv = cstwtp+freight+insuran+communic+techsupp+accnting+inhouse+salescom+othserv
		replace serv=. if serv<=0
		replace rserv=. if rserv<=0

		gen GO=output/10
		gen RGO=routput/10
		replace RGO=. if RGO<=0

		gen VA=valadded/10
		gen RVA=rva/10
		drop rva
		replace RVA=. if RVA<=0

		gen WI=(energ+matls+serv)/10
		gen RI=(renerg+realmats+rserv)/10
		gen RN=renerg/10
		gen RV=rserv/10
		gen RM=realmats/10
		gen RMV=RM+RV
		gen RNV=RN+RV
		gen WN=energ/10
		gen WV=serv/10
		gen WM=matls/10
		gen WMV=WM+WV
		gen WNV=WN+WV

		gen L = totlab
		gen WL = totwages/10
		gen WLS= skwages/10
		gen WLU = unskwages/10


		gen cpi1985=cpi/16.6436
		gen RL=WL/cpi85
		gen RLS=WLS/cpi85
		gen RLU=WLU/cpi85

		gen K = rcapstock/10
		replace K=. if K<=0

		gen netinv = rinvest/10
		replace netinv=. if netinv<=0
		gen nominv=invnbld + invnmach + invnveh + invubld + invumach + invuveh - salesbld -salesuma - salesuve + impbld + impmach + impveh + prodbld + prodmach + prodveh
		replace nominv=nominv/10
		replace nominv=. if nominv<=0

		gen si = log(WI/GO)
		gen sl = log(WL/GO)
		gen sn = log(WN/GO)
		gen snv = log(WNV/GO)

		replace defoutput=defoutput/100

		//=================================================================================================================================================
		//=================================================================================================================================================
		// Specific to each industry
		//=================================================================================================================================================
		//=================================================================================================================================================
		//replace si=. if si>1.5
		//replace K=. if K>2500
		//=================================================================================================================================================
		//=================================================================================================================================================
		//=================================================================================================================================================
		//=================================================================================================================================================

		gen export=.
		replace export=0 if exports~=. & year>1989
		replace export=1 if exports~=. & exports>0
		gen import=.
		replace import=0 if rawmatsi~=. & rawmatsi==0
		replace import=1 if rawmatsi~=. & rawmatsi>0

		gen wages=.
		replace wages=totwages/L
		bys year: egen mm=median(wages) if wages~=.
		gen hiwag=.
		replace hiwag=0 if wages<mm & wages!=.
		replace hiwag=1 if wages>mm & wages!=.
		drop mm

		keep id year GO RGO VA RVA RN WN RV WV RM WM RMV WMV RNV WNV RI WI RL WL RLS WLS RLU WLU K netinv nominv L si sl sn snv defoutput export import hiwag adv ratiowcbc ///
				faland stassets defcons_ye defcons_ye XX152 XX147 XX148 XX149 XX150 XX151 netprof

		reshape wide GO RGO VA RVA RN WN RV WV RM WM RMV WMV RNV WNV RI WI RL WL RLS WLS RLU WLU K netinv nominv L si sl sn snv defoutput export import hiwag adv ratiowcbc faland stassets defcons_ye XX152 XX147 XX148 XX149 XX150 XX151 netprof, i(id) j(year)
		reshape long GO RGO VA RVA RN WN RV WV RM WM RMV WMV RNV WNV RI WI RL WL RLS WLS RLU WLU K netinv nominv L si sl sn snv defoutput export import hiwag adv ratiowcbc faland stassets defcons_ye XX152 XX147 XX148 XX149 XX150 XX151 netprof, i(id) j(year)

		bys year: egen def=mean(defoutput)
		drop defoutput

		rename netinv INV

		foreach i of varlist GO RGO VA RVA RN RV RM RMV RNV RI RL RLS RLU K INV nominv L def{
			gen log`i'=log(`i')
		}	

		compress
		sort id year
		save miss_data_jan_28_2011, replace

		bys year: egen tot=sum(RGO)
		gen s=RGO/tot
		gen temp=s*ln(RGO)
		bys year: egen q_index=sum(temp)
		drop temp
		gen temp=s*RGO
		bys year: egen q_index_l=sum(temp)
		drop tot s temp

		foreach i of varlist GO RGO VA RVA RN RV RM RMV RNV RI RL RLS RLU K INV nominv L def{
			replace log`i'=-1000.0 if `i'==. | `i'<=0
		}
		 
		foreach i of varlist GO RGO VA RVA RN WN RV WV RM WM RMV WMV RNV WNV RI WI RL WL RLS WLS RLU WLU K INV nominv L {
			gen I`i'=1
			replace I`i'=0 if `i'==. | `i'<=0
			replace `i'=-1000.0 if I`i'==0
		}	

		foreach i of varlist si sl sn snv {
			gen I`i'=1
			replace I`i'=0 if `i'==.
			replace `i'=-1000.0 if I`i'==0
		}

		replace export=-1000 if export==.
		replace import=-1000 if import==.
		bys year: egen mm=median(logK) if logK>-900
		gen sizek=.
		replace sizek=0 if logK<mm & logK>-999
		replace sizek=1 if logK>=mm
		drop mm
		replace sizek=-1000 if sizek==.

		bys year: egen mm=median(ratiowcbc) 
		gen wcbc=.
		replace wcbc=0 if ratiowcbc<=mm
		replace wcbc=1 if ratiowcbc>mm
		drop mm ratiowcbc
		replace wcbc=-1000 if wcbc==.

		gen l_q_index=ln(q_index)
		gen l_q_index_l=ln(q_index_l)
		drop def logdef q_index q_index_l

		replace hiwag=-1000 if hiwag==.
		replace adv=-1000 if adv==.

		#delimit;
		order GO logGO IGO RGO logRGO IRGO VA logVA IVA RVA logRVA IRVA si Isi sl Isl sn Isn snv Isnv RN logRN IRN WN IWN RV logRV IRV WV IWV 
			RM logRM IRM WM IWM RMV logRMV IRMV WMV IWMV RNV logRNV IRNV WNV IWNV
			RI logRI IRI WI IWI RL logRL IRL WL IWL 
			RLS logRLS IRLS WLS IWLS RLU logRLU IRLU WLU IWLU K logK IK id year 
			INV logINV IINV nominv lognominv Inominv L logL IL l_q_index l_q_index_l 
			export import sizek hiwag adv wcbc; 
		#delimit cr
		sort id year
		compress
		
		//This command is my own; I don't know why log(0) was set to -1000, but I undo it here;
		//since this has to be done to match GNR's OLS results, I'm pretty sure they do the same thing
		recode logL logK logRI logRGO (-1000 = .)

end









use "dta/ChileAnalysis2.dta", clear
formdata `ind'

gen time = year
gen yg = logRGO
gen k = logK
gen l = logL
gen i = logRI

recode K L RI (-1000 = .)

replace si = . if si == -1000

gen ll=l*l
gen kk=k*k
gen ii=i*i
gen lk=l*k
gen li=l*i
gen ki=k*i
gen lki=l*k*i
*gen lll=ll*l
*gen kkk=kk*k
*gen iii=ii*i


/************************************************************/
/**  The code below runs the share regression using a      **/ 
/**  a (log)polynomial approximation.  This step     	   **/
/**  recovers the output elasticity of flexible inputs 	   **/
/**                                                        **/
/**  The initial values are set using an OLS regression    **/
/**  of shares (si) on a polynomial in inputs.  If the     **/
/**  non-linear least squares procedure (nl) fails to      **/
/**  converge, the initial values can be changed.	   **/
/************************************************************/

// Obtain initial values from OLS refression of flexible inputs shares on all productive inputs //
regress si l k i ll kk ii lk li ki lki if si~=. & l~=. & k~=. & i~=.
matrix test = e(b)
predict temp_c if si~=. & l~=. & k~=. & i~=.
replace temp_c = temp_c - _b[_cons]

//Correct constant term for initial values //
egen mtemp_c = min(temp_c)
scalar ntemp_c=mtemp_c
drop temp_c mtemp_c
scalar ntemp_c=-ntemp_c + 0.1

// Using initial values from OLS, regress shares using log polynomial approximation //
#delimit;
nl ( si = ln({g0=ntemp_c} + {gl=_b[l]}*l + {gk=_b[k]}*k + {gi=_b[i]}*i + {gll=_b[ll]}*ll + 
   {glk=_b[lk]}*lk + {gli=_b[li]}*li + {gkk=_b[kk]}*kk + {gki=_b[ki]}*ki +
   {gii=_b[ii]}*ii + {glki=_b[lki]}*lki) ) if si~=. & l~=. & k~=. & i~=., iter(100);
#delimit cr

// Construct predicted elasticities of output with respect to intermediate inputs //
predict ielas if l~=. & k~=. & si~=. & i~=.

// Perform mean adjustment of predicted elasticities //
predict eg if l~=. & k~=. & si~=. & i~=., resid
replace eg=-eg
egen mexp_eg=mean(exp(eg))
replace ielas=ielas-ln(mexp_eg)
replace ielas=exp(ielas)

// Rename and save mean adjusted coefficients of estimation //
mat beta=e(b)
svmat double beta
ren beta1 g0
ren beta2 gl
ren beta3 gk
ren beta4 gi
ren beta5 gll
ren beta6 glk
ren beta7 gli
ren beta8 gkk
ren beta9 gki
ren beta10 gii
ren beta11 glki
foreach var of varlist g0-glki {
	egen s`var'=mean(`var')
	drop `var'
	ren s`var' `var'
	replace `var' = `var' / mexp_eg
}	
clear matrix

// Compute the integral //
gen integ_G_I = g0+gl*l+gk*k+gll*ll+gkk*kk+glk*lk + (gi*i+gli*li+gki*ki+glki*lki)/2 + gii*ii/3
replace integ_G_I=integ_G_I*i 

// Construct script_Y //
gen vg = yg - eg - integ_G_I

// Generate lagged variables //
tset id time
gen vg_1=L.vg
gen l_1=L.l
gen k_1=L.k
gen ll_1=L.ll
gen kk_1=L.kk
gen lk_1=L.lk
gen vg_2=L2.vg
gen l_2=L2.l
gen k_2=L2.k
gen ll_2=L2.ll
gen kk_2=L2.kk
gen lk_2=L2.lk

/************************************************************/
/**  Now to recover the remaining coefficients associated  **/
/**  with capital and labor (the constant of the PDE)      **/
/**                                                        **/
/**  The initial values are set using an OLS regression    **/
/**  of vg on log capital and loglabor                     **/
/************************************************************/

// Set initial values using OLS estimation of script_Y on capital, labor, their squares and interaction term //
reg vg l k ll kk lk 
matrix test2 = e(b)
matrix test3 = test2[1,1..5]
matrix test3[1,1] = test2[1,"l"]
matrix test3[1,2] = test2[1,"k"]
matrix test3[1,3] = test2[1,"ll"]
matrix test3[1,4] = test2[1,"kk"]
matrix test3[1,5] = test2[1,"lk"]
matrix drop test2

// Recover capital and labor coefficients using gmm procedure //
// gmm_prod function computes moment conditions using productivity innovation and capital and labor inputs //
gmm gmm_prod if vg!=. & l!=. & k!=. & vg_1!=. & l_1!=. & k_1!=., one nequations(5) parameters(al ak all akk alk) from(test3) winitial(identity) rhs(vg l k ll kk lk vg_1 l_1 k_1 ll_1 kk_1 lk_1) conv_maxiter(100)


mat beta=e(b)
svmat double beta
ren beta1 al
ren beta2 ak
ren beta3 all
ren beta4 akk
ren beta5 alk
foreach var of varlist al-alk {
	egen s`var'=mean(`var')
	drop `var'
	ren s`var' `var'
}	

/************************************************************/
/**  Generate productivity in logs and levels, and compute **/
/**  the output elasticities of labor and capital.         **/
/************************************************************/

// generate log productivity //
gen logomega=vg-al*l-ak*k-all*ll-akk*kk-alk*lk
gen omega=exp(logomega)

// Generate output elasticities of labor and capital inputs //
gen lelas=gl*i+2*gll*li+glk*ki + gli*ii/2 + glki*ki*i/2 + al + 2*all*l + alk*k
gen kelas=gk*i+2*gkk*ki+glk*li + gki*ii/2 + glki*li*i/2 + ak + 2*akk*k + alk*l

sum *elas






/*

. recode K L RI (-1000 = .)

. sum mexp_eg

    Variable |       Obs        Mean    Std. Dev.       Min        Max
-------------+--------------------------------------------------------
     mexp_eg |     52272     1.03742           0    1.03742    1.03742

. gen prod = 1.03742* omega
(32898 missing values generated)

. reg K l(1/2).K l(1/2).L l(1/2).prod

      Source |       SS       df       MS              Number of obs =   14243
-------------+------------------------------           F(  6, 14236) =       .
       Model |  3.7963e+13     6  6.3272e+12           Prob > F      =  0.0000
    Residual |  7.7170e+11 14236  54207563.3           R-squared     =  0.9801
-------------+------------------------------           Adj R-squared =  0.9801
       Total |  3.8735e+13 14242  2.7198e+09           Root MSE      =  7362.6

------------------------------------------------------------------------------
           K |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
           K |
         L1. |   1.597243    .007808   204.57   0.000     1.581938    1.612547
         L2. |  -.6069227   .0079927   -75.93   0.000    -.6225896   -.5912559
             |
           L |
         L1. |   10.94184   1.042077    10.50   0.000      8.89923    12.98444
         L2. |  -1.851277   1.085344    -1.71   0.088    -3.978692    .2761382
             |
        prod |
         L1. |   .0226414   .5581861     0.04   0.968    -1.071476    1.116759
         L2. |   1.363307   .5392621     2.53   0.011     .3062832    2.420331
             |
       _cons |  -2079.842   422.8294    -4.92   0.000    -2908.642   -1251.041
------------------------------------------------------------------------------

. reg L l(1/2).K l(1/2).L l(1/2).prod

      Source |       SS       df       MS              Number of obs =   14232
-------------+------------------------------           F(  6, 14225) =15561.09
       Model |   337573596     6  56262266.1           Prob > F      =  0.0000
    Residual |  51431529.2 14225  3615.57323           R-squared     =  0.8678
-------------+------------------------------           Adj R-squared =  0.8677
       Total |   389005126 14231   27335.052           Root MSE      =   60.13

------------------------------------------------------------------------------
           L |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
           K |
         L1. |   .0006951   .0000638    10.90   0.000     .0005701    .0008201
         L2. |  -.0005811   .0000653    -8.90   0.000    -.0007091   -.0004532
             |
           L |
         L1. |   .7863375    .008508    92.42   0.000     .7696608    .8030143
         L2. |   .1478766    .008858    16.69   0.000     .1305138    .1652394
             |
        prod |
         L1. |   .0310601   .0045602     6.81   0.000     .0221215    .0399986
         L2. |   -.010386   .0044043    -2.36   0.018    -.0190189    -.001753
             |
       _cons |  -19.91258   3.449735    -5.77   0.000    -26.67451   -13.15064
------------------------------------------------------------------------------

. 



gen logs = log(stassets)

. reg logs l(1/2).logs

      Source |       SS       df       MS              Number of obs =    1625
-------------+------------------------------           F(  2,  1622) = 1339.01
       Model |  5431.48625     2  2715.74312           Prob > F      =  0.0000
    Residual |  3289.69712  1622  2.02817332           R-squared     =  0.6228
-------------+------------------------------           Adj R-squared =  0.6223
       Total |  8721.18337  1624   5.3701868           Root MSE      =  1.4241

------------------------------------------------------------------------------
        logs |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
        logs |
         L1. |   .4930425   .0249598    19.75   0.000     .4440856    .5419994
         L2. |   .3094285   .0247868    12.48   0.000      .260811     .358046
             |
       _cons |   1.983175   .1196155    16.58   0.000     1.748558    2.217792
------------------------------------------------------------------------------


*/
